### Alizade, Dancygier, Ditlmann 
### "National Penalties Reversed"
### Replication Code 
### Figure A5
### For questions, contact jalizade@princeton.edu

# empty environment
rm(list = ls())


#setwd("")
setwd("C:/Users/Jey/Dropbox/WZB/NaturalizationExperiment/Submission/JOP/replication_JOP/data")


# load necessary packages
library(readstata13)
library(ggplot2)
library(margins)

# load data set
dat <- read.dta13("data_experimental.dta")


### regressions ###

# center-right
right_mod <- glm(e1_response ~ e1_treat_turkish + e1_treat_dual + e1_treat_benefits, dat[dat$e1_leftpol==0,], family = "binomial")

# center-left
left_mod <- glm(e1_response ~ e1_treat_turkish + e1_treat_dual + e1_treat_benefits, dat[dat$e1_leftpol==1,], family = "binomial")


### AMEs ###

# center-right
ame_right <- margins(right_mod)

# center-left
ame_left <- margins(left_mod)


### Figure ###

# create data frame with estimates
ests <- rbind.data.frame(summary(ame_right)[,c(1:3)], 
                         summary(ame_left)[,c(1:3)])
ests$factor <- rep(c("Benefits", "Dual", "Turkish"), 2)
ests$pol <- rep(0:1, each=3)
ests$pol <- factor(ests$pol, levels = c(1,0))

# plot
fig1 <- ggplot(ests, aes(x = factor, y = AME), group = pol)
fig1 <- fig1 + geom_hline(yintercept = 0, linetype = "dashed")
fig1 <- fig1 + geom_errorbar(aes(ymax = AME + SE, ymin = AME - SE,
                                 group = pol, color = factor(pol)), size = 1.2, 
                             position = position_dodge(width = 0.5), width = 0)
fig1 <- fig1 + geom_errorbar(aes(ymax = AME + SE*1.96, ymin = AME - SE*1.96,
                                 group = pol, color = factor(pol)), size = 0.4, 
                             position = position_dodge(width = 0.5), width = 0.2)
fig1 <- fig1 + geom_point(aes(x = factor, y = AME, group = pol, shape = factor(pol), fill = factor(pol), color = factor(pol)), 
                          position = position_dodge(width = 0.5), size = 3)
fig1 <- fig1 + xlab("") + ylab("Average Marginal Effect")
fig1 <- fig1 + coord_flip()
fig1 <- fig1 + scale_shape_manual(values = c(24, 22), name = "Politician\nPartisanship", labels=c("Center-Left", "Center-Right"))
fig1 <- fig1 + scale_fill_manual(values = c("red", "black"), name = "Politician\nPartisanship", labels=c("Center-Left", "Center-Right"))
fig1 <- fig1 + scale_color_manual(values = c("red", "black"), name = "Politician\nPartisanship", labels=c("Center-Left", "Center-Right"))
fig1 <- fig1 + theme(legend.title=element_text(size=10, face = "bold"))
fig1 <- fig1 + guides(shape = guide_legend(reverse = T), color = guide_legend(reverse = T), fill = guide_legend(reverse = T))
fig1
